Exploratory Data Analysis

A Demonstration of Tools

Michael Clark https://m-clark.github.io
2020-07-08

Table of Contents


Introduction

In R there are many tools available to help you dive in and explore your data. However, in consulting I still see a lot of people using base R’s table and summary functions, followed by a lot of work to get the result into a more presentable format. My own frustrations led to me creating a package (tidyext) for personal use in this area. While that suits me fine, there are tools that can go much further with little effort. Recently, Staniak & Biecek Staniak and Biecek (2019) wrote an article in the R Journal exploring several of such packages, so I thought I’d try them out for myself, and take others along with me for that ride.

Note that these are first impressions, and I haven’t really dived deeply into any of the packages, and may be missing some key features (apologies to the package authors!). But that is also part of the point for this sort of thing. These aren’t modeling packages, and we have a good idea of what we want in EDA, so these should be easy to pick up and use.

Packages

I’ve updated the article’s table 1, which adds roughly a year’s worth of downloads.

package downloads debut
janitor 495548 2016-10-03
summarytools 199235 2014-08-11
DataExplorer 186456 2016-03-01
visdat 184922 2017-07-11
funModeling 102122 2016-02-07
arsenal 76190 2016-12-30
dlookr 44482 2018-04-27
dataMaid 42799 2017-01-02
inspectdf 23351 2019-04-24
xray 19409 2017-11-22
RtutoR 16272 2016-03-12
ExPanDaR 15464 2018-05-11
exploreR 14459 2016-02-10
SmartEDA 13674 2018-04-06
explore 10125 2019-05-16

A more informative assessment of usage would be in average monthly downloads.

package average_monthly_downloads
janitor 10772.783
visdat 4997.892
DataExplorer 3518.038
summarytools 2767.153
funModeling 1891.148
arsenal 1731.591
dlookr 1588.643
inspectdf 1459.438
dataMaid 995.326
explore 675.000
xray 588.152
ExPanDaR 572.741
SmartEDA 488.357
RtutoR 307.019
exploreR 267.759

Here is a visualization of their growth over time.

Selected packages

Arbitrary criteria

I’ll outline my reasons for selecting some packages to explore and not others. These reasons are somewhat, but not necessarily, arbitrary, and may leave out some viable newer packages. For the data scenario, I am assuming messy data of the sort that might have hundreds of columns of mixed data types, potentially with lots of missingness, attributes that are only applicable to subsets of the data (e.g. branching logic in surveys), etc.1

Here is my criteria for selection:

Things I’m not as concerned about:

In the end, I would like a package that is well put together and will make common tasks easier for me and potentially save me time in creating reports/presentation.

Conceptual organization

Staniak & Biecek note two general phases of data exploration, each with specific tasks, based on the CRISP-DM standard Wirth and Hipp (2000).2

In the first place, I will focus on tools for understanding, particularly description and validity, as they refer to exploration tasks solely as visualization, which is a perk, but something I’m more inclined to do myself.

I’m definitely less interested in the preparation though I will speak about it more towards the end. ‘Cleaning’ in the article refers to mean/median imputation, something I’ve never bothered to do for reasons that have been noted in the statistical literature for a very long time. The other transformations are easy, and probably should be more explicitly documented in your code. Furthermore, in creating derived attributes, things like category merging and standardization depend on the subset of the data used, so it would probably be better to be more explicit than automatic. Also, if things like an automated PCA is viable for your situation, it probably is very simple data (i.e. all variables of the same type), in which case, most of these tools probably won’t add much value to you anyway.

So with that in mind here are the ones I will explore (in alphabetical order):

The Data

I’ve chosen the heart disease data originally available from the UCI repository (object name hd). It contains a mixture of data types but isn’t too unwieldy, as it’s already been cleaned and has few columns (it’s from my noiris package). For our purposes, I’ve additionally added some random missingness, and created an hd_sample which has only a couple columns to cut down on the display.

Data Description

For data description, we are interested in things like, dimensions of the data, variable types, and maybe even meta meta-data like the object’s size in RAM. I also add univariate data summaries as it’s something you’d always want to know and usually report, though these are thought of as data exploration in the article’s delineation.

Preliminary

To give a sense of what my preferences are, consider my own functions. The following actually calls a separate numerical summary function as well as a categorical variable function, and both return ‘tidy’ data frames that can immediately be used for presentation (e.g. via kableExtra) and visualization (e.g. ggplot2), or drilling down to only selections of the output.


library(tidyext)

describe_all(hd, extra = T)  # also describe_all_num, describe_all_cat

$`Numeric Variables`
# A tibble: 8 x 12
  Variable          N   Mean    SD   Min    Q1 Median    Q3   Max `% Missing` Distinct Zeros
  <chr>         <dbl>  <dbl> <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl>       <dbl>    <dbl> <dbl>
1 age             286  54.3   8.97    29  47.2   55    61    77             6       40     0
2 resting_bp      286 132.   17.7     94 120    130   140   200             6       49     0
3 cholesterol     290 246.   52.3    126 211.   240.  274   564             4      148     0
4 resting_ecg     286   0.52  0.52     0   0      1     1     2             6        3   141
5 max_heartrate   292 150.   22.8     71 136.   153   167.  202             4       90     0
6 old_peak        289   1.03  1.16     0   0      0.8   1.6   6.2           5       38    95
7 n_vessels       287   0.72  1.02     0   0      0     1     4             5        5   167
8 heart_disease   284   0.53  0.5      0   0      1     1     1             6        2   133

$`Categorical Variables`
# A tibble: 22 x 4
   Variable            Group            Frequency   `%`
   <chr>               <fct>                <int> <dbl>
 1 sex                 male                   199   66 
 2 sex                 female                  92   30 
 3 sex                 <NA>                    12    4 
 4 chest_pain_type     typical angina         135   45 
 5 chest_pain_type     non-anginal pain        84   28.
 6 chest_pain_type     atypical angina         47   16 
 7 chest_pain_type     asymptomatic            21    7.
 8 chest_pain_type     <NA>                    16    5 
 9 fasting_blood_sugar lt_120                 242   80 
10 fasting_blood_sugar gt_120                  44   15 
# … with 12 more rows

I also have options, such as the following.


hd %>% 
  select(age, sex, heart_disease) %>% 
  describe_all(
  digits = 2,
  include_NAcat = FALSE,   # don't include NA as a category
  include_numeric = TRUE,  # allows numeric variables with few levels given in max_levels argument as categorical
  max_levels = 3,
  sort_by_freq = TRUE,
  extra = TRUE
)

$`Numeric Variables`
# A tibble: 2 x 12
  Variable          N  Mean    SD   Min    Q1 Median    Q3   Max `% Missing` Distinct Zeros
  <chr>         <dbl> <dbl> <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl>       <dbl>    <dbl> <dbl>
1 age             286 54.3   8.97    29  47.2     55    61    77           6       40     0
2 heart_disease   284  0.53  0.5      0   0        1     1     1           6        2   133

$`Categorical Variables`
# A tibble: 4 x 4
  Variable      Group  Frequency   `%`
  <chr>         <fct>      <int> <dbl>
1 sex           male         199    68
2 sex           female        92    32
3 heart_disease 1            151    53
4 heart_disease 0            133    47

For grouped output, I can use the underlying functions. Here we look at summaries for a couple numerical variables by sex.


hd %>% 
  num_by(
    main_var  = vars(age, cholesterol),
    group_var = sex,
    extra = TRUE
  )

# A tibble: 6 x 13
# Groups:   sex [3]
  sex    Variable        N  Mean    SD   Min    Q1 Median    Q3   Max `% Missing` Distinct Zeros
  <chr>  <chr>       <dbl> <dbl> <dbl> <dbl> <dbl>  <dbl> <dbl> <dbl>       <dbl>    <dbl> <dbl>
1 male   age           189  53.6   8.9    29  47      54   59      77           5       38     0
2 male   cholesterol   188 239.   43     126 210.    234. 268.    353           6      112     0
3 female age            86  55.3   9.2    34  49.2    56   62      76           7       32     0
4 female cholesterol    90 262.   64.7   141 215.    253  300.    564           2       75     0
5 <NA>   age            11  58.3   7.9    44  54      61   63.5    68           8        9     0
6 <NA>   cholesterol    12 245.   60.4   166 210     237  261.    394           0       12     0

The categorical variable functionality is very similar.


hd %>% 
  cat_by(
    main_var  = chest_pain_type,
    group_var = sex
  )

# A tibble: 13 x 5
# Groups:   sex [3]
   sex    chest_pain_type      N `% of Total` `% of sex`
   <chr>  <chr>            <int>        <dbl>      <dbl>
 1 female asymptomatic         3        0.990       3.26
 2 female atypical angina     14        4.62       15.2 
 3 female non-anginal pain    31       10.2        33.7 
 4 female typical angina      34       11.2        37.0 
 5 female <NA>                10        3.30       10.9 
 6 male   asymptomatic        18        5.94        9.05
 7 male   atypical angina     32       10.6        16.1 
 8 male   non-anginal pain    49       16.2        24.6 
 9 male   typical angina      94       31.0        47.2 
10 male   <NA>                 6        1.98        3.02
11 <NA>   atypical angina      1        0.330       8.33
12 <NA>   non-anginal pain     4        1.32       33.3 
13 <NA>   typical angina       7        2.31       58.3 

As these are tidy tibbles, they are essentially ready for presentation.


describe_all_num(hd) %>% 
  kableExtra::kable()
Variable N Mean SD Min Q1 Median Q3 Max % Missing
age 286 54.26 8.97 29 47.25 55.0 61.00 77.0 6
resting_bp 286 131.85 17.66 94 120.00 130.0 140.00 200.0 6
cholesterol 290 246.32 52.27 126 211.25 240.5 274.00 564.0 4
resting_ecg 286 0.52 0.52 0 0.00 1.0 1.00 2.0 6
max_heartrate 292 149.83 22.82 71 135.50 153.0 167.25 202.0 4
old_peak 289 1.03 1.16 0 0.00 0.8 1.60 6.2 5
n_vessels 287 0.72 1.02 0 0.00 0.0 1.00 4.0 5
heart_disease 284 0.53 0.50 0 0.00 1.0 1.00 1.0 6

These functions serve most of my needs for initial peeking at the data. They return a tibble/data.frame class object that makes for easy presentation and visualization. The underlying code uses the widely used and tested tidyverse packages, and mostly adheres to standard programming conventions. This is the same sort of thing I’m looking for.

arsenal

We begin alphabetically with the arsenal package. Here we use tableby for a generic summary as well as grouped summary. The result is markdown, so for presentation in an R Markdown document, one must use the chunk option results='asis'.


library(arsenal)
hd_sample %>%
  tableby( ~ ., data = .) %>%
  summary()
Overall (N=303)
age
   N-Miss 17
   Mean (SD) 54.262 (8.973)
   Range 29.000 - 77.000
cholesterol
   N-Miss 13
   Mean (SD) 246.317 (52.271)
   Range 126.000 - 564.000
resting_bp
   N-Miss 17
   Mean (SD) 131.850 (17.663)
   Range 94.000 - 200.000
sex
   N-Miss 12
   female 92 (31.6%)
   male 199 (68.4%)
chest_pain_type
   N-Miss 16
   asymptomatic 21 (7.3%)
   atypical angina 47 (16.4%)
   non-anginal pain 84 (29.3%)
   typical angina 135 (47.0%)
heart_disease
   N-Miss 19
   Mean (SD) 0.532 (0.500)
   Range 0.000 - 1.000

hd_sample %>%
  tableby(sex ~ ., data = .) %>%
  summary()
female (N=92) male (N=199) Total (N=291) p value
age 0.143
   N-Miss 6 10 16
   Mean (SD) 55.279 (9.187) 53.566 (8.871) 54.102 (8.989)
   Range 34.000 - 76.000 29.000 - 77.000 29.000 - 77.000
cholesterol < 0.001
   N-Miss 2 11 13
   Mean (SD) 261.456 (64.737) 239.144 (43.028) 246.367 (52.017)
   Range 141.000 - 564.000 126.000 - 353.000 126.000 - 564.000
resting_bp 0.284
   N-Miss 4 13 17
   Mean (SD) 133.682 (19.732) 131.199 (16.934) 131.996 (17.882)
   Range 94.000 - 200.000 94.000 - 192.000 94.000 - 200.000
chest_pain_type 0.106
   N-Miss 10 6 16
   asymptomatic 3 (3.7%) 18 (9.3%) 21 (7.6%)
   atypical angina 14 (17.1%) 32 (16.6%) 46 (16.7%)
   non-anginal pain 31 (37.8%) 49 (25.4%) 80 (29.1%)
   typical angina 34 (41.5%) 94 (48.7%) 128 (46.5%)
heart_disease < 0.001
   N-Miss 2 16 18
   Mean (SD) 0.744 (0.439) 0.426 (0.496) 0.531 (0.500)
   Range 0.000 - 1.000 0.000 - 1.000 0.000 - 1.000

Here is an example of categorical-only output using freqlist.


with(hd, table(sex, chest_pain_type)) %>% 
  freqlist() %>% 
  summary()
sex chest_pain_type Freq Cumulative Freq Percent Cumulative Percent
female asymptomatic 3 3 1.09 1.09
atypical angina 14 17 5.09 6.18
non-anginal pain 31 48 11.27 17.45
typical angina 34 82 12.36 29.82
male asymptomatic 18 100 6.55 36.36
atypical angina 32 132 11.64 48.00
non-anginal pain 49 181 17.82 65.82
typical angina 94 275 34.18 100.00

Other options include a function for doing pairwise comparisons for a repeated measure3, comparing datasets, and doing a bunch of bivariate regressions.

Pros

The tableby summary is essentially a ‘Table 1’, which is an unfortunately named display of descriptive stats for a sample of a given study4. My clients often want that, and I’ve actually used a package specifically geared towards providing it just to save me the headache (tableOne), so I definitely find that aspect useful. However, Table 1’s almost invariably have needless statistical output or analysis, and are not a model for reporting that I would choose to go for. As you can see, the layout is not going to be viable for more than a few variables, though it is common practice to present them that way in journal articles regardless of verbosity/legibility.

Issues

I’m not thrilled with markdown output, as I have little control over it, but it’s fine. Also, the default layout is not so succinct, and it appears it’s trying to emulate SAS functions, which are not models for presentation in my opinion. One can use a function to write the output of a single table to html, but I’d rather it just be amenable to the document I’m already creating, or create a report for me of all tables. Lastly, I’m not crazy about using the summary function to get the output. The underlying list objects are not ‘print ready’, so using summary as an argument with default of TRUE would make more sense to me design-wise.

DataExplorer

Now we move to DataExplorer. Let’s introduce the package with the introduce function. I already like this, as it provides good info that extends what you’d get with str or glimpse.


library(DataExplorer)

introduce(hd)

# A tibble: 1 x 9
   rows columns discrete_columns continuous_colu… all_missing_col… total_missing_v… complete_rows total_observati…
  <int>   <int>            <int>            <int>            <int>            <int>         <int>            <int>
1   303      14                6                8                0              214           150             4242
# … with 1 more variable: memory_usage <dbl>

It also can display this information visually in two different ways. I like this, but can’t say I’d ever have a reason to actually use it.


plot_intro(hd)


plot_str(list(hd = hd, gapminder = gapminder_2019))

We can focus on the missingness, which is nice, but I don’t find a use for bar simple bar charts, as they actually make what is only a single row of information harder to parse, and without any visual interest.


plot_missing(hd_sample)

DataExplorer can also plot distributions, e.g. bar and histograms, for the actual values of the categorical/discrete variables.


plot_bar(hd)

But this sort of thing only takes a couple lines of ggplot to do on your own, which you can then customize far more easily. At least DataExplorer can save you the trouble and sorts the result in an ordered fashion by default.


hd %>% 
  select_if(is.character) %>% 
  pivot_longer(everything(), names_to = 'variable', values_to = 'value') %>% 
  ggplot(aes(x = value)) +
  geom_bar() +
  coord_flip() +
  facet_wrap(~ variable, scales = 'free') 

Similarly there are QQ plots, scatterplots and more. I always try to get people to try correlation plots in lieu of large correlation matrices, and DataExplorer provides this. The nice thing is that it will automatically create indicator variables for levels of categorical variables, but beyond that there are issues. For one, it’s diagonal is reversed from the typical presentation of correlation matrices5. If you don’t drop the missing, the plot isn’t as useful, because it will include NA as a factor level, but to its credit DataExplorer has an option to only focus on continuous or discrete output. In addition, ‘centered’ is a rather odd choice for alignment of the x axis in my opinion, so this would take additional work to make presentable.


plot_correlation(
  hd,
  ggtheme = theme_minimal(),
  cor_args = list("use" = "pairwise.complete.obs")
)

There are packages that do this specifically, such as heatmaply. Before these came around I already had my own function, that also does an internal factor analysis to sort the variables, and produces an interactive result. So while so this aspect of DataExplorer might be useful to others, it doesn’t appeal much to me, especially given the other issues.


hd %>% 
  select_if(is.numeric) %>% 
  cor(use = 'pair') %>% 
  visibly::corr_heat()

If the various options of output, or only some of them, appeal to you, they can all be nicely wrapped up in an automatic report. Various settings for each type of output can be set with an additional function (configure_report) or just passing a list of arguments, including which functions to use, ggplot2 theme, etc. You can see the report here, but I show a screenshot.


create_report(
  hd,
  y = 'heart_disease',
  output_dir   = 'other_docs',
  output_file  = 'data_explorer_report.html',
  report_title = 'My Data Description'
)

Pros

I think many would like at least some functionality in DataExplorer, as well as many of the visualizations. The ease with which to generate a report should also be sufficient for anyone’s personal use, and with some tweaking, presentation to others. It also uses data.table under the hood, so likely can handle large data with efficiency.

Not covered here, but DataExplorer also has functionality for feature processing and engineering, for example, collapsing sparse categories, dummy coding, etc.

Issues

The issues I have are pretty minor with this one aside from unnecessary statistical analysis and visualization choices. I would also recommend pryr::object_size rather than base R’s function.

SmartEDA

I wanted to look at the SmartEDA package because the figure in the article was of a clean report. Let’s start our exploration with the basic ExpData function.


library(SmartEDA)
ExpData(hd)

                                Descriptions       Obs
1                         Sample size (Nrow)       303
2                    No. of Variables (Ncol)        14
3                   No. of Numeric Variables         8
4                    No. of Factor Variables         0
5                      No. of Text Variables         6
6                   No. of Logical Variables         0
7                    No. of Unique Variables         0
8                      No. of Date Variables         0
9   No. of Zero variance Variables (Uniform)         0
10     %. of Variables having complete cases    0% (0)
11 %. of Variables having <50% missing cases 100% (14)
12 %. of Variables having >50% missing cases    0% (0)
13 %. of Variables having >90% missing cases    0% (0)

We can look at ExpNumStat to get some basic stats for numeric variables.


ExpNumStat(hd, round = 1)

          Vname Group  TN nNeg nZero nPos NegInf PosInf NA_Value Per_of_Missing     sum min   max  mean median   SD  CV
1           age   All 303    0     0  286      0      0       17            5.6 15519.0  29  77.0  54.3   55.0  9.0 0.2
3   cholesterol   All 303    0     0  290      0      0       13            4.3 71432.0 126 564.0 246.3  240.5 52.3 0.2
4 max_heartrate   All 303    0     0  292      0      0       11            3.6 43751.0  71 202.0 149.8  153.0 22.8 0.2
5      old_peak   All 303    0    95  194      0      0       14            4.6   296.9   0   6.2   1.0    0.8  1.2 1.1
2    resting_bp   All 303    0     0  286      0      0       17            5.6 37709.0  94 200.0 131.8  130.0 17.7 0.1
   IQR Skewness Kurtosis
1 13.8     -0.2     -0.5
3 62.8      1.2      4.5
4 31.8     -0.6      0.0
5  1.6      1.3      1.7
2 20.0      0.7      0.9

We can look at ExpNumStat to get some basic stats for grouped output also, and set various options.


ExpNumStat(
  hd_sample,
  by = "GA",
  gp = "sex",
  Qnt = c(.1, .9),
  Outlier = TRUE,
  round = 1
)

         Vname      Group  TN nNeg nZero nPos NegInf PosInf NA_Value Per_of_Missing   sum min  max  mean median   SD
1          age    sex:All 303    0     0  286      0      0       17            5.6 15519  29   77  54.3   55.0  9.0
4          age   sex:male 199    0     0  189      0      0       10            5.0 10124  29   77  53.6   54.0  8.9
7          age sex:female  92    0     0   86      0      0        6            6.5  4754  34   76  55.3   56.0  9.2
10         age     sex:NA   0    0     0    0      0      0        0            NaN     0 Inf -Inf   NaN     NA   NA
2  cholesterol    sex:All 303    0     0  290      0      0       13            4.3 71432 126  564 246.3  240.5 52.3
5  cholesterol   sex:male 199    0     0  188      0      0       11            5.5 44959 126  353 239.1  234.5 43.0
8  cholesterol sex:female  92    0     0   90      0      0        2            2.2 23531 141  564 261.5  253.0 64.7
11 cholesterol     sex:NA   0    0     0    0      0      0        0            NaN     0 Inf -Inf   NaN     NA   NA
3   resting_bp    sex:All 303    0     0  286      0      0       17            5.6 37709  94  200 131.8  130.0 17.7
6   resting_bp   sex:male 199    0     0  186      0      0       13            6.5 24403  94  192 131.2  130.0 16.9
9   resting_bp sex:female  92    0     0   88      0      0        4            4.3 11764  94  200 133.7  132.0 19.7
12  resting_bp     sex:NA   0    0     0    0      0      0        0            NaN     0 Inf -Inf   NaN     NA   NA
    CV  IQR Skewness Kurtosis   10% 90% LB.25% UB.75% nOutliers
1  0.2 13.8     -0.2     -0.5  42.0  66   26.6   81.6         0
4  0.2 12.0     -0.2     -0.5  42.0  65   29.0   77.0         0
7  0.2 12.8     -0.2     -0.7  42.0  66   30.1   81.1         0
10  NA   NA      NaN      NaN    NA  NA     NA     NA         0
2  0.2 62.8      1.2      4.5 186.9 309  117.1  368.1         5
5  0.2 57.2      0.1     -0.3 184.7 299  124.4  353.4         0
8  0.2 85.0      1.4      4.2 197.0 331   87.8  427.8         1
11  NA   NA      NaN      NaN    NA  NA     NA     NA         0
3  0.1 20.0      0.7      0.9 110.0 152   90.0  170.0         9
6  0.1 20.0      0.7      0.6 110.0 152   90.0  170.0         4
9  0.1 20.0      0.7      0.9 108.0 160   90.0  170.0         5
12  NA   NA      NaN      NaN    NA  NA     NA     NA         0

There is also some visualization of relationships with a given target variable.


ExpNumViz(hd_sample, target = 'cholesterol')

[[1]]


[[2]]

Here are the default categorical data summaries. This is a nice and clean data.frame presentation.


ExpCTable(hd_sample) 

          Variable            Valid Frequency Percent CumPercent
1              sex           female        92   30.36      30.36
2              sex             male       199   65.68      96.04
3              sex               NA        12    3.96     100.00
4              sex            TOTAL       303      NA         NA
5  chest_pain_type     asymptomatic        21    6.93       6.93
6  chest_pain_type  atypical angina        47   15.51      22.44
7  chest_pain_type               NA        16    5.28      27.72
8  chest_pain_type non-anginal pain        84   27.72      55.44
9  chest_pain_type   typical angina       135   44.55      99.99
10 chest_pain_type            TOTAL       303      NA         NA
11   heart_disease                0       133   43.89      43.89
12   heart_disease                1       151   49.83      93.72
13   heart_disease               NA        19    6.27      99.99
14   heart_disease            TOTAL       303      NA         NA

As there was with numeric variables, there is also visualization for the categorical variables.


ExpCatViz(hd_sample) 

[[1]]


[[2]]


[[3]]

As far as reporting, SmartEDA also provides this functionality. There are some options you can fiddle with, like using your own template, changing the default theme, etc.


ExpReport(
  hd,
  theme   = visibly::theme_clean(),
  op_dir  = 'other_docs/',
  op_file = 'smarteda.html'
)

Pros

SmartEDA is fairly intuitive to use. It returns a data frame, and some of the less verbose output is quite ready to go. It can also generate a report in automatic fashion.

Issues

A lot of this isn’t very useful to me, such as parallel coordinate plots, ‘outlier’ analysis, etc. I’m not crazy about the naming conventions, both for the functions and arguments (every single function in the package begins with Exp). The default color schemes for some output is simply not viable for presentation, and the report doesn’t really have any options for controlling the output like DataExplorer did. I altos had issues with attempting to get it to work outside of the current working directory for the initial document creation. And changing the sn or sc options, which I’m not sure why I would only want a sample of plots rather than specifying which plots I’d want specifically, didn’t appear to actually change the resulting document.

summarytools

The summarytools package provides four main functions to work with, but I’m going to skip those and go straight to the tool that uses those key functions and puts their results into a very nice presentation.


library(summarytools)

dfSummary(
  hd,
  varnumbers = FALSE,
  round.digits = 2,
  plain.ascii = FALSE,
  style = "grid",
  graph.magnif = .33,
  valid.col = FALSE,
  tmp.img.dir = "img"
)

Data Frame Summary

hd

Dimensions: 303 x 14
Duplicates: 0

Variable Stats / Values Freqs (% of Valid) Graph Missing
age
[numeric]
Mean (sd) : 54.3 (9)
min < med < max:
29 < 55 < 77
IQR (CV) : 13.8 (0.2)
40 distinct values 17
(5.61%)
sex
[character]
1. female
2. male
92 (31.6%)
199 (68.4%)
12
(3.96%)
chest_pain_type
[character]
1. asymptomatic
2. atypical angina
3. non-anginal pain
4. typical angina
21 ( 7.3%)
47 (16.4%)
84 (29.3%)
135 (47.0%)
16
(5.28%)
resting_bp
[numeric]
Mean (sd) : 131.8 (17.7)
min < med < max:
94 < 130 < 200
IQR (CV) : 20 (0.1)
49 distinct values 17
(5.61%)
cholesterol
[numeric]
Mean (sd) : 246.3 (52.3)
min < med < max:
126 < 240.5 < 564
IQR (CV) : 62.8 (0.2)
148 distinct values 13
(4.29%)
fasting_blood_sugar
[character]
1. gt_120
2. lt_120
44 (15.4%)
242 (84.6%)
17
(5.61%)
resting_ecg
[numeric]
Mean (sd) : 0.5 (0.5)
min < med < max:
0 < 1 < 2
IQR (CV) : 1 (1)
0 : 141 (49.3%)
1 : 142 (49.6%)
2 : 3 ( 1.1%)
17
(5.61%)
max_heartrate
[numeric]
Mean (sd) : 149.8 (22.8)
min < med < max:
71 < 153 < 202
IQR (CV) : 31.8 (0.2)
90 distinct values 11
(3.63%)
exer_angina
[character]
1. no
2. yes
190 (66.2%)
97 (33.8%)
16
(5.28%)
old_peak
[numeric]
Mean (sd) : 1 (1.2)
min < med < max:
0 < 0.8 < 6.2
IQR (CV) : 1.6 (1.1)
38 distinct values 14
(4.62%)
slope
[character]
1. flat
2. negative
3. positive
133 (45.6%)
138 (47.3%)
21 ( 7.2%)
11
(3.63%)
n_vessels
[numeric]
Mean (sd) : 0.7 (1)
min < med < max:
0 < 0 < 4
IQR (CV) : 1 (1.4)
0 : 167 (58.2%)
1 : 60 (20.9%)
2 : 37 (12.9%)
3 : 18 ( 6.3%)
4 : 5 ( 1.7%)
16
(5.28%)
defect
[character]
1. fixed_defect
2. normal
3. reversible_defect
160 (56.1%)
17 ( 6.0%)
108 (37.9%)
18
(5.94%)
heart_disease
[numeric]
Min : 0
Mean : 0.5
Max : 1
0 : 133 (46.8%)
1 : 151 (53.2%)
19
(6.27%)

This is exactly what I want- basic, not overwhelming and redundant information, a usable data frame object, simple visuals to enhance the output without adding to the total information you have to parse (and which can be turned off), and basic categorical information. In one function, I have pretty much all I’d need, but with control to tweak as necessary.


library(summarytools)

hd_sample %>% 
  group_by(sex) %>% 
  dfSummary(
    varnumbers = FALSE,
    round.digits = 2,
    plain.ascii = FALSE,
    na.col = FALSE,
    style = "grid",
    graph.magnif = .33,
    valid.col = FALSE,
    tmp.img.dir = "/tmp"
  )

Data Frame Summary

hd_sample

Group: sex = female
Dimensions: 92 x 6
Duplicates: 0

Variable Stats / Values Freqs (% of Valid) Graph
age
[numeric]
Mean (sd) : 55.3 (9.2)
min < med < max:
34 < 56 < 76
IQR (CV) : 12.8 (0.2)
32 distinct values
cholesterol
[numeric]
Mean (sd) : 261.5 (64.7)
min < med < max:
141 < 253 < 564
IQR (CV) : 85 (0.2)
75 distinct values
resting_bp
[numeric]
Mean (sd) : 133.7 (19.7)
min < med < max:
94 < 132 < 200
IQR (CV) : 20 (0.1)
33 distinct values
sex
[character]
1. female 92 (100.0%)
chest_pain_type
[character]
1. asymptomatic
2. atypical angina
3. non-anginal pain
4. typical angina
3 ( 3.7%)
14 (17.1%)
31 (37.8%)
34 (41.5%)
heart_disease
[numeric]
Min : 0
Mean : 0.7
Max : 1
0 : 23 (25.6%)
1 : 67 (74.4%)

Group: sex = male
Dimensions: 199 x 6
Duplicates: 0

Variable Stats / Values Freqs (% of Valid) Graph
age
[numeric]
Mean (sd) : 53.6 (8.9)
min < med < max:
29 < 54 < 77
IQR (CV) : 12 (0.2)
38 distinct values
cholesterol
[numeric]
Mean (sd) : 239.1 (43)
min < med < max:
126 < 234.5 < 353
IQR (CV) : 57.2 (0.2)
112 distinct values
resting_bp
[numeric]
Mean (sd) : 131.2 (16.9)
min < med < max:
94 < 130 < 192
IQR (CV) : 20 (0.1)
43 distinct values
sex
[character]
1. male 199 (100.0%)
chest_pain_type
[character]
1. asymptomatic
2. atypical angina
3. non-anginal pain
4. typical angina
18 ( 9.3%)
32 (16.6%)
49 (25.4%)
94 (48.7%)
heart_disease
[numeric]
Min : 0
Mean : 0.4
Max : 1
0 : 105 (57.4%)
1 : 78 (42.6%)

Group: sex = NA
Dimensions: 12 x 6
Duplicates: 0

Variable Stats / Values Freqs (% of Valid) Graph
age
[numeric]
Mean (sd) : 58.3 (7.9)
min < med < max:
44 < 61 < 68
IQR (CV) : 9.5 (0.1)
44 : 1 ( 9.1%)
46 : 1 ( 9.1%)
54 : 2 (18.2%)
60 : 1 ( 9.1%)
61 : 2 (18.2%)
62 : 1 ( 9.1%)
65 : 1 ( 9.1%)
66 : 1 ( 9.1%)
68 : 1 ( 9.1%)
cholesterol
[numeric]
Mean (sd) : 245.2 (60.4)
min < med < max:
166 < 237 < 394
IQR (CV) : 50.8 (0.2)
12 distinct values
resting_bp
[numeric]
Mean (sd) : 128.5 (11.6)
min < med < max:
110 < 127.5 < 142
IQR (CV) : 20 (0.1)
110 : 1 ( 8.3%)
117 : 1 ( 8.3%)
120 : 4 (33.3%)
135 : 1 ( 8.3%)
138 : 1 ( 8.3%)
140 : 3 (25.0%)
142 : 1 ( 8.3%)
sex
[character]
All NA’s
chest_pain_type
[character]
1. atypical angina
2. non-anginal pain
3. typical angina
1 ( 8.3%)
4 (33.3%)
7 (58.3%)
heart_disease
[numeric]
Min : 0
Mean : 0.5
Max : 1
0 : 5 (45.5%)
1 : 6 (54.5%)

Pros

The summarytools package provides a single function that produces a ready-to-present table within whatever document I’m already using. Very nice!

Issues

The only nitpicky stuff I have with this package is that the ctable function isn’t useful to me, and some of the numeric description is not really necessary. I had issues with the trying to use the suggested tmp directory for the image files, but found it easier to just use a project folder.

Data Validity

As a reminder, data validity is more about providing checks on the data rather than summarizing it per se. These are packages more geared toward spotting or dealing with data issues.

dataMaid

The primary utility of dataMaid is examination of data consistency, but it will also provide basic summaries too. In the end it creates a useful report as well. I came across this via the author’s presentation at the Ann Arbor R User Group, and have actually used it before for a project with success. The result can be found here, though here is a screenshot.


library(dataMaid)

dataMaid::makeDataReport(
  hd,
  output = 'html',
  file = 'other_docs/dataMaid_report',
  replace = TRUE,
  maxDecimals = 1
) 

As you can see, it provides an overall summary of the data frame including variable types, unique values, missing percentage, and potential problems. Problems mostly regard outliers, which is fine to inspect, but arbitrarily set.

Pros

dataMaid is highly customizable, but the default is already a great way to get a good sense of whether your data is doing what it’s supposed to do.

Issues

I initially couldn’t get anything but a pdf document even though I clearly specified html as per the documentation, so I’m not sure what’s going on there, but I eventually sorted it out. The default outlier check, while better than most would do, seemed too sensitive, but this is nitpicky. In the past I’ve found it be slow in rendering the output for a larger data set, but this isn’t something you would usually have to run more than once after getting your options squared away.

janitor

Sadly, I still regularly receive data from Excel (and SPSS), which means I have to worry about things that I shouldn’t have to worry about, like whether dates will actually be treated appropriately, if the column names are usable, or whether the imported file contains 10000 empty rows because someone accidentally hit the space bar.

The janitor package provides a few simple tools that many will probably find useful at some point, and I regularly use its remove_empty function after importing anything from Excel. As an example I’ll add an empty column with an Excel-like name, a duplicated row, and a useless column with a constant value (all very typical in Excel files I receive).


library(janitor)

hd2 = hd

hd2$`test this - ridiculous (name)` = NA           # a terribly formatted name that has no values
hd2$`why.is-this_here` = 'same value everywhere'   # constant column
hd2[nrow(hd2) + 1,] = hd2[1,]                      # add a duplicated row

hd2 %>% 
  clean_names() %>% 
  colnames()

 [1] "age"                       "sex"                       "chest_pain_type"           "resting_bp"               
 [5] "cholesterol"               "fasting_blood_sugar"       "resting_ecg"               "max_heartrate"            
 [9] "exer_angina"               "old_peak"                  "slope"                     "n_vessels"                
[13] "defect"                    "heart_disease"             "test_this_ridiculous_name" "why_is_this_here"         

And here is the remove_* functionality.


hd3 = hd2 %>% 
  remove_empty() %>% 
  remove_constant() 

# useless columns/rows removed 
colnames(hd3)

 [1] "age"                 "sex"                 "chest_pain_type"     "resting_bp"          "cholesterol"        
 [6] "fasting_blood_sugar" "resting_ecg"         "max_heartrate"       "exer_angina"         "old_peak"           
[11] "slope"               "n_vessels"           "defect"              "heart_disease"      

c(nrow(hd2), nrow(hd3))

[1] 304 304

Just a couple little things like that are very useful. Otherwise, there is functionality for dates, tables, etc. that some might find useful in a pinch, but one would probably will have more with the other summary functions in packages we’ve demoed, and more advanced users may simply just use lubridate, stringr and similar packages directly.

Pros

Provides useful functionality not seen in the other packages.

Issues

Not a whole lot going on here relative to some of the other packages, but what is there is likely to be useful to many.

visdat

The package visdat6 is, as its name implies, purely for visualization, and this includes missingness, correlation, and more. We can start by visualizing the missing data.


library(visdat)
vis_dat(hd_sample)

We can look at variable types along with the missingness. I like the integer/double mix, which might point to an issue if all the data should be integer valued.


vis_guess(hd)

We can visualize the correlations as we did with DataExplorer, and default is pairwise. This is one of the cleaner correlation plots I’ve come across, but there is no way to order it meaningfully.


vis_cor(hd %>% select_if(is.numeric))

One nice feature I haven’t seen elsewhere is to visualize a given expected value across a data set.


hd %>% 
  select(sex) %>% 
  vis_expect(expectation =  ~ .x == 'male')

We can also compare whole data sets. Here the differences are with missing values in one and not the other7.


hd2 = noiris::heart_disease %>% 
  select(colnames(hd_sample)) %>% 
  mutate_if(is.factor, as.character)

vis_compare(
  hd_sample,
  hd2
)

Pros

This package adds some useful functionality we haven’t seen. In addition, the underlying code adheres to open science standards. There is a complementary package naniar for more ways to deal with missing data.

Issues

This isn’t for numeric description, so can only supplement typical EDA as we’ve had before. You may have to do some pre-processing for some functions unlike with other tools (e.g. subsetting to numeric). But these are minor issues at best given the package’s purpose.

Data Cleaning and Transformations

As mentioned previously, I’m not interested in using these packages for doing this, and based on Staniak and Biecek, only a couple of them provide this, DataExplorer being the only package that we’ve explored here. Mean/median/mode imputation is a great way to attenuate correlation, so I’m not interested in doing that. I don’t even know what an ‘outlier’ is outside of a modeling framework, so I’m definitely not interested in doing something about extreme values based on a univariate analysis8. Discretizing numeric variables should almost never be done9, so that functionality is not desirable. And between base R functions like scale and log, or the rescale function from scales, lumping factors and similar via forcats, I don’t really need another package for this sort of thing.

Summary

Staniak & Biecek Staniak and Biecek (2019) provide a nice summary of many packages that can be used for automated exploratory data analysis. I wanted to explore and demonstrate some of these in more detail. My impression is that many of these packages have notably useful functionality, though a lot of it may be overlapping, or superfluous from a more serious analytical standpoint. In addition, I personally wouldn’t use the word automated to describe most of the functionality for most of these packages, as multiple specific functions, possibly with different options for each, would need to be used to ultimately generate a data driven product. That’s not necessarily a bad thing, as one should know what’s in the data if you want to make a decision based on it.

Here is a rough list of features I made based on quick inspections.

Table 1: See Table 4 in Staniak & Biecek for more details.
Package Ready to use output Code Quality Visualization Report Generation Dedicated Website
arsenal 👍
DataExplorer 👍 👍 👍 👍 👍
SmartEDA 👍 👍 👍 👍
summarytools 👍 👍
visdat 👍 👍 👍 👍
dataMaid 👍 👍
janitor 👍 👍 👍

In general I would imagine I’d use summarytools or DataExplorer, with visdat and janitor to fill in additional information or help with the data processing. Interestingly, these also appear to be the most popular packages under consideration10.

Description packages

arsenal

DataExplorer

SmartEDA

summarytools

Validation packages

visdat

dataMaid

janitor

Exercises

Find a data set you like and do the following. I suggest using a manageable one or some subset so that you don’t get an overwhelming amount of output.

Staniak, Mateusz, and Przemyslaw Biecek. 2019. “The Landscape of R Packages for Automated Exploratory Data Analysis.” The R Journal. https://journal.r-project.org/archive/2019/RJ-2019-033/index.html.

Wirth, Rüdiger, and Jochen Hipp. 2000. “CRISP-Dm: Towards a Standard Process Model for Data Mining.” In Proceedings of the 4th International Conference on the Practical Applications of Knowledge Discovery and Data Mining, 29–39. Springer-Verlag London, UK.


  1. However, for our demonstration I’ll be using something a little more wieldy.↩︎

  2. To be honest, I wasn’t familiar with the CRoss Industry Standard Process for Data Mining until reading the article citing it. I don’t get the impression any particular methodology is actually consciously thought about by the vast majority of practicing data scientists, but it’s useful in providing a framework for the content here.↩︎

  3. The data we’re using isn’t paired and their documentation doesn’t provide a working example.↩︎

  4. I see ‘Table 1’ used mostly by folks in medical fields, who subsequently place it as table 2, 3 or whatever as would normally be the case.↩︎

  5. I suspect this may have more to do with ggplot than what is desired. But I also found that the orientation/alignment of names was different for the report versus the inline presentation.↩︎

  6. The article by Staniak & Biecek places visdat with description rather than validity, but the data frame comparison, expectation investigation, and missing exploration seem more validity issues than data summaries.↩︎

  7. The arsenal package also provided this functionality. It’s quite verbose, so I didn’t show it, but if you need a more in depth comparison it can be recommended.↩︎

  8. Outliers are an indicator of model weakness/limitation/failure, not a data problem.↩︎

  9. The only times I do this is with an already discrete variable, for example, going from 7 values to 5 values, or coarsening that might be applied in very large data situations for computational efficiency. Issues with discretizing are noted here, and while there are better ways to do it if you feel you must, you probably won’t have seen it done in practice, and the chosen method likely wouldn’t hold across repeated data.↩︎

  10. The janitor package, which was not explored in the original article, is actually notably more popular than any of the other packages.↩︎